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The influence of dark matter (DM) on the growth of supermassive black holes 
(SMBHs) is studied. It is shown that gravitational scattering of DM particles on bulge 
stars leads to diffusion of DM in phase space {m, m^,/} (m denotes the angular mo- 
mentum and I is the radial action). Appropriate diffusion coefHcients are calculated for 
different bulge models, and it is argued that the diffusion along m axis is the most impor- 
tant effect. It is shown that this process leads to noticeable flow of DM into the black hole 
(BH), resulting in its power-law growth: Mi,h oc t®/^^. Comparison with observational 
data shows that, in principle, this effect may explain observed masses of SMBHs. Special 
attention is paid to the corrections related to the innermost region of BH gravitational 
influence and the diffusion along I axis. Their influence on the BH growth law is shown 
to be negligible. 



1. Introduction 

The interaction of a supermassive black hole in a galaxy center with dark matter halo 
is already much investigated -"^ ^ . There an adiabatic invariant approach is applied 
for different halo structure. The initial halo profiles are taken to be self-similar 
(power-law) p ^ 7-^"(0 < a < 2), isothermal (p ~ r~^) and NFW profileJ^! [p ~ 
—I — T-rh-/ — w)- In Ref. ilithe method of adiabatic invariant is used to calculate small 
changes in orbital parameters of particles caused by slow variation of gravitational 
potential due to the black hole growth. This approach has two drawbacks: absorption 
by BH is not taken into account, and appropriate change of DM distribution function 
is neglected. The latter means that the loss cone is always full, thus leading to 
overestimation of dark matter flow. 

More correct approach is used in Refs. El El There the evolution of DM distri- 
bution function due to absorption by BH and changes in BH mass and loss cone 
parameters are considered consistently. The authors come to the conclusion that 
for current values of black hole masses the fraction of dark matter inside them is 
rather small. 
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This approach was developed in Ref. where the change of dark matter distri- 
bution function outside the loss cone due to diffusion in phase space was taken into 
account. It was demonstrated that this diffusion caused by gravitational scatter- 
ing of DM particles on stars effectively refills the loss cone and determines the BH 
growth law. This effect was shown to give reasonable estimate for observed black 
holes masses. 

The present work further develops the mentioned approach by accounting for 
three other factors: change of gravitational potential in the vicinity of black hole, 
modification of star distribution in this region, and the conditions under which the 
diffusion becomes effectively two-dimensional. 

The paper is organized as follows. In the first section we introduce our model 
of dark matter halo and the galactic bulge. In the second section we write down 
the kinetic equation which describes the evolution of dark matter halo in phase- 
space, and calculate several quantities related to star distribution and dark matter 
motion. Calculation of the diffusion coefficients is performed in the third section. 
Also estimated are the diffusion timescales and corrections related to the central 
region. In the fourth section the solution of the diffusion equation in one-dimensional 
form is presented, the black hole growth law is established, and the corrections to the 
DM flow due to two-dimensional diffusion and the central region are investigated. It 
is shown that these corrections in most cases play no significant role in the diffusion 
process. Finally, in the section 5 we present some numerical results and discuss 
them. 

1.1. Dark matter in galaxies 

While the most part of matter in the Universe is now proved to be dark, the nature 
of the dark matter still remains unclear. In most models, however, the significant 
or dominant part of dark matter is taken to consist of cold dark matter (CDM), i.e. 
non-relativistic particles interacting only gravitationally. This leads to the absence 
of thermodynamical equilibrium in DM gas and the necessity of kinetic approach. 
Following this approach it was shown that the initial inhomogenities in DM distribu- 
tion grow and form spherically-symmetric structures - non-dissipative gravitational 
singularities (NGSs) of different scales, having the same internal structure . Gen eral 
analytic theory of NGS formation 1^ (consistent with numerical simulations 1 ^ 1 ^ ^ 1 and 
other analytic theories E|) predicts that all these structures have similar density 
profiles described by the following formula: 

pir)^ Ar-^ , C = 12/7. (1) 

The density profiles of DM haloes appearing in numerical simulations (see Refs. 0] 
have different forms, typically shallower in the centre. While being more realistic 
in the outer parts than the profile Q due to hierarchical formation mechanism, they 
cannot describe well the very inner region becau se of l imited spatial resolution, 
thereby allowing different values of the inner slope ElMI ranging from —1 to —1.5, 
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all of which are within the errors i n th e simulations. Furthermore, observations do 
not exclude cuspy density profiles li^, but still are unable to determine the inner 
slope precisely. Therefore, the adopted profile does not contradict observations and 
is convenient for a number of reasons (simplicity of distribution function and tight 
relation between angular momentum and energy). 

As initial inhomogenities are not spherically symmetric, the particles of NGS 
possess angular momenta (though total angular momentum of NGS is zero) . High 
degree of contraction at the nonlinear stage makes these momenta relatively small 
compared to maximal possible angular momenta for corresponding orbital sizes. It 
turns out that angular momentum of a particle is related to its radial action / by 

the formula: = Iq^P, where / = ^ Vr{r)dr, Vr = y^2(£;- *(r)) - ^ is 
the radial velocity, ^'(r) is the gravitational potential, r± are the turnpoints, and 
— 0.1 is a small parameter. 
In this case, the distribution function of DM in canonical variables /, m, has 
the following form: 

/(/,m,m,) = /o/i/'5(m2-/o'/') (2) 

The average value of m is zero, though the average is nonzero. 

The formation of NGS of galactic scale is followed by the formation of galaxy 
itself due to infall of baryonic matter into gravitational well of NGS, its subsequent 
cooling and formation of galactic structures - disc, bulge and halo. After this stage 
the total gravitational potential is no longer defined by DM; conversely, in the 
central region it is mainly directed by baryons. Galaxy formation process is slow 
enough as compared to the dynamical time of DM particles, so that the gravitational 
potential evolution is adiabatic. This fact ensures that radial action is integral of 
motion (adiabatic invariant). It is defined by equation 



1 r+ / m2 
I{E,m)^- drp{E-^{r))-— (3) 

where r_, r+ are minimal and maximal distances of a DM particle from center; 
r_ <C due to small angular momentum (m <^ I). 

Our main task is to study the central region of a galaxy - the bulge, which we 
may consider to be spherically symmetric. Under this assumption angular momen- 
tum is also conserved, and we obtain that the DM distribution function written in 
variables /, m, does not depend on time. Spatial distribution of DM is given by 
relation 

p(r) (2^)3— — / dm dm, dE^ . .f(-E,m), (4) 



47rr^ 



E-'i'-^ 



is the energy interval where the expression under the radical sign is non-negative. 
Let us assume that the density and the potential of the bulge have power-law 
profiles: 

n{r) cx r-(2"") , Mb{r) cx r"+i , *(r) cx r" , 
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Then one can show that the radial action may be approximated (with the 
accuracy better than 8%) in a factorised form 



I = Io{r+)C{p), M 



m 



(5) 



1/2 



where function C(/i) at small argument /i is linear: 

« C(0) - , C(0)>0, (6) 
Hence the DM density profile is also power-law: 

p(r) = A'r-«', e' = -f + ^a (7) 

Let us now overview the bulge structure and proceed to estimation of DM frac- 
tion in the bulge. 



1.2. Bulge overview 

Bulge is the central p art of a galaxy, for our Galaxy its radius being about 1 kpc and 
mass about IQ^'^Mq ^1 We assume it to be spherically symmetric. Bulge density 
profile can be derived from the dependency of stars velocity dispersion on distance 
to galactic center; we take it to be power-law: n{r) — jybr"^^""^. For isothermal 
bulge, i.e. having uniform velocity dispersion, a = 0. Milky Way bulge may be 
approximated with velocity dispersion a oc r^/"* for r < 50 pcl^^, which corresponds 
a = 0.5. 

It should be noticed that currently available observations of distant galaxies 
usually do not have e noug h resolution to derive velocity dispersion on scales less 
than a hundred parsecs'-'^. We assume their bulges to be isothermal in central parts, 
but one can show that the exact form of density profile does not matter much; only 
its integral characteristics affect the final result. 

Generally, in the center of a galaxy a compact object is located, which is assumed 
to be a supermassive black hole (for a recent review see Ref. I17|l . For our own 
Galaxy it is now proved that the compact object Sgr A* is a black hole with mass 
h'hh — 2.9 • IO^A/q . The innermost star cluster has the following density profile: 



n(r) = no 



-3/2 



, where ng — 10® 



pc'^ 



, ro = 0.02 pc 



(8) 



M{r) = /([ Af^n(r') 47rr'2 dr' = 
^MgHoro'^ ( ) • ^'^^ ^^'^ sake of simplicity we have taken the masses of all 



Hence the mass of this cluster 

3/2 
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stars to be equal Mg = Mq. Let us define the influence radius Rq of BH by the 
condition that the total mass of stars inside this radius equals that of the BH: 

/ j^j 2 \ 2 / 3 

Ro = I ^TT I w 1 pc at the moment. (9) 

So we may assume that inside this radius the gravitational potential is dominated 
by BH, while outside it the influence of BH is negligible. 

Notice that BH and bulge evolution changes values of no and Rq. If we assume 
that the star density is continuous at Rq while Rq increases due to BH growth, then 
the following relation between Albh and Rq arises: 

A4/^ = YMs»7b^o'+". (10) 

Here fjb is the coefficient in the bulge density profile; for isothermal bulge with 
velocity dispersion ctq fjt ~ 2tt'g m (^"^^ ^^)- 

1.3. Dark matter in the bulge 

The total mass of Milky Way's dark matter halo can be assessed to be Mh ^ 
IO^^Mq, its radius to be Rh ~ 100 kpcEH. From © we derive adiabatic invariant 
for particles with _R+ = Rh to be Imax = O.S5{G Mh RhY^'^ ■ If we suppose that 
the initial distribution can be extrapolated up to i?+ = Rh, then tha total halo 



and hence 



Mh = ^i27:ffoI,na.'^' , (11) 
1 

^ ^-^(2^ (Gi?^jW • (^2) 

Then the value of /o for our Galaxy is 6 • I0^g{^)-^/^. 

In this case the relative amount of DM inside the bulge is given by 

^ _ MpM _ Mh (G Mb Rb f'^ _( Mh \ ^ / \ ^ , . 

Mb {GMHRnf"' Mb \Mb ) \Rh ) ^ ' ^ ' 

This means that for the model distribution Q the mass of dark matter inside bulge 
is comparable to the baryonic mass. The same relation is valid in the innermost 
region (BH domain of influence). However, the distribution Q is applicable only 
in some vicinity of the galactic center, and hence the actual values of /o and T can 
deviate from the assessments (|T^ IT^ . 

Observations and modelling of mass distribution in the Galaxy do not give tight 
restrictions on the amount of DM in the bulge. Depending on the model, Tf may 
comprise from 1/4 to 1/3 of the total mass inside the central 3 kpc regional but 
not the major fraction of the bulge mass. So, for the final calculations of the BH 
mass we should decrease the value of /□ H12|l by a factor of 3-4 or more. 
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2. Problem definition and particle motion parameters 

Our goal is to determine the possible fraction of DM in BH, and to derive the growth 
law of BH due to absorption of DM. Firstly, one can prove that the direct capture of 
particles with momenta less than nig = 'iGAhh joes not lead to significant growth 
of BH since their mass calculated from distribution |(2J) with current value of is 
several orders of magnitude less than Mth 

The simplest process which changes the distribution function of DM is gravita- 
tional scattering of DM particles on stars identical to Coulomb scattering in plasma. 



2.1. Kinetic equation 

The dark matter distribution function satisfies the following equation 
df{r, V, t) 



dt 



{H^J} = St{f} 



(14) 



where Hq is the Hamilton function for gravitational interaction, {,} is the Poisson 
bracket, St{f} is the collision term in Landau form 



st{f} 



d{^iVi) 



Bi 



J' 



d{fj,Vj 



B^j cPv' 



(15) 



Here entities with primes refer to stars and other to DM particles; masses are 
denoted as ^, /x' to avoid confusion with momentum, u — v — if is the relative 
velocity, Lc ^ 15 is the Coulomb logarithm. Since the mass of a star Mg = fi' ^ iJ,, 
we neglect the first term in equation: 



st{f} 



d_ 



dvi 



(16) 



f'iv',r)d\' 



As soon as collision frequency of DM particles with stars is much less than their 
orbital motion frequency, one can rewrite the equation (|14|l in the form averaged 
over a period in action-angle variables {Ik,4>k}', {Ik} = {I Tm^mz}'. 

df{{Ik},t) 



dt 



St{f} 



Then the equation 1)16(1 for collision term looks as follows: 



st{f} 



d 



Rkl 



Rkl 



(27r)3 



df 

dll_ 

dh dh 



OVi OVj 



(17) 

(18) 
(19) 



The star distribution around a SMBH and their diffusion towards the disruption 
boundary was studied in R,efs. l2,'^ll24l Two-dimensional Fokker-Planck equation was 
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written in coordinates E and m, and solved analytically with certain simplifications, 
introducing the notion of energy-dependent loss-cone. A fraction of stars with small 
angular momenta is eliminated during one dynamical time, and the diffusion along 
energy axis leads to capture of stars with critical (minimal) energy. We use radial 
action / instead of energy, because diffusion along / axis does not lead to capture; 
therefore the absorption process looks more clear. Also, our distribution function is 
time-dependent and the black hole mass changes in time. 

Observations show that the distribution function (DF) of stars may be consid- 
ered as isotropic, i.e. not depending on momentum, even in the central region 
For isotropic DF the tensor Wij takes the following formal 



ViV 



= A{E, r),5y - B{E, r) ^ , (20) 

OO 

A ^ ^-^GHi-L.^ I dE'f'iE') { ; ' ; : (21) 




*(r) 
00 

B = '-^GHI^^L^ I dE'nE') {],,,^lz t (22) 



r-)\E>E' 

From 10 © it follows that adiabatic invariant may be represented as follows: 

I{E,m)f=i.J{E)-hc,m. (23) 

Note that in the case of Coulomb (a = —1) and oscillator (a = 2) potential the 
constant ba — 1, and in the case of isothermal potential (or close to it) ba ~ 0.6. 

It is convenient to perform linear variable change: {I,m,mz} — > {J,m,mz} 
Since it is linear, the expression (|19|l for tensor Rki does not change. In addition, 
due to small value of parameter Iq, the initial distribution function (0) has the same 
form in new variables: 

MJ,m) = foJ'/'Sim^-lo'j^) (24) 

The initial distribution function does not depend on m^, hence the solution of 
kinetic equation H17(l will not depend on neither. We rewrite the expression (|18() 
for collision term as follows: 

St{f} = —^r—m i?22^ V ^12^ + TTT M?12t; ^ -^11^ ' (25) 



mom \ dm oJ J oJ \ dm a J 

where the diffusion coefficients (|19|l are the following^ 

R..= [%)\{A-By) , (26) 
i?i2- (H) {{A-B)m) , 
i?22 ^[at^^B- 
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Here the averaging over angle variables is defined as 




T = 2r-, «. = A/2(i?-*(r))-^, (28) 

where T is the particle oscillation period, Vj. is its radial velocity. 

Before proceeding to the solution of the kinetic equation ifTTI 1^ one should 
define the distribution function of stars. 



2.2. Distribution of stars 

We assume that the density profile of stars in the bulge is power-law: n{r) = 

(nl^) ' ^'^'^ ^^'^ gravitational potential is also power-law ^'(r) = ^o'^"- In this 
case isotropic di stri bution function of stars can be written as a power-law depen- 
dence on energy 

fiv' ,t) ^ F^E-f" , i?=^ + *(r), (29) 

n{r) = / /'(«', r) dh'' = / Fo(-- + ^o^")"'' 47rw'2 dv' = 
Jo Jo 2 



= il(27r)3/2FoVl/oi-''r"(i-« , hence 7 = (/?- (30) 

r(/3) 2 

We consider two particular cases: a) self-consistent potential of stars in the bulge 
(from observations '^'^ it follows that the star distribution in the bulge is close to 
isothermal, so we take the power-law index 7 close to 2) and b) the central region 
of bulge where potential is dominated by the black hole. 

In the first case we have 

— — = = GMsUorf^r^ ^ , hence 

dr r'^ — 7 



a = (Tor"/^ is the velocity dispersion, = ^0 2(i-a) i (32) 
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In the limit a — > these expressions describe isothermal star cluster. Let us 
rewrite them in this particular case separately: 

v|/(r) = In — , ^'o = 2(70^ , cq is the velocity dispersion, (34) 

f'{E) = Fo exp(— ^) , Fo = no(2wo')"'/' ■ (36) 

In the second case the particle dynamics is governed by the Coulomb potential 
of the black hole. 

= ; (37) 

r 

Restricting our consideration to the case 7 = 3/2 which follows from observations 
11^, we obtain from (|30|l /3 = 0. Hence the distribution function does not depend on 
energy and equals Fq at E < 0: 

where Mth and Rq are related by l(Tn|l . 



2.3. Dark matter motion parameters 

In this section we present the motion parameters for DM particles for practically 
important cases. 



1 /"'+ / m2 
I(E, m) — — \ 2(E — ^'(r)) — dr - adiabatic invariant; 

TT V 



I{E, 0) = J{E) , r± - turnpoints; 

T{E,m) — 27r— particle oscillation period; 

oE 

r_ m 

= Xmin < 1 , — - = /i < 1 . 

r+ J 
For power-law potential (|31|) : 

i+f 1 



V"' 

T(i?) = i±l ^^^^ , (40) 
1 — a (To 



2fTo 1 — a / y 47r(l — a) 



E a \ ' _ I a 
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In the limit a ^ we obtain formulas for isothermal potential H34|l : 
1 



J{E) 



-.(jf)rj, 



T{E) - V^r+/ao , 
/ E 



ro exp 



Xmin 



For Coulomb potential ll^ : 



T{E)^ 



TT GM, 



bh 



V2 \EV/^ ' 
GMbh 



\E\ 



Xmin — ^ 



(42) 
(43) 
(44) 

(45) 
(46) 
(47) 



3. The diffusion coefficients 

Two opposite cases should be distinguished: for particles with apocenter distances 
r+(J) < Rq {Rq is the radius of BH's domain of influence we calculate coef- 

ficients using Coulomb potential, and for other particles we use the expressions for 
self-consistent bulge potential (this corresponds to assumption that these particles 
spend the most part of orbital period outside BH's domain of influence). 

3.1. The diffusion coefficients for bulge 

Below we present the quantities A, A — B for power-law (including isothermal) star 
distribution: 



A 



A-B 



3\/27r 
4 



3V27r 



where 

Acii) 



(?)V2^ + l)-('3-i) {)2 



P-1 

aia + 1)^ r(/3) 



2q I 2'^ ' 2' 2q 



r(/3-|) ' 



(«72'?+l) 



(3-1 



-f('- ii-l--"^ 

2 2' 2^ 



a (To"' 



a{r) ' 



m ' 2' 2q 



a(a + l)V? r(/?) 



(48) 
(49) 



2 r(/3-|)' 

is the hypergeometric function. 
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and for isothermal profile 



My) 



exp(- 



exp(- 



Hx) = Jo exp{-*-j)dt 
is the probability integral. 



The graphs of functions Aa{x), Q,a{x) are displayed on fig^ 
Now let us calculate the diffusion coefficients J?i2, R22 according to 
First of all, we are interested in low values of momentum, so we neglect the 
second term in R22, which is ~ fi^ times smaller than the first one. 



R2 



3V2tt 



GMsL,r+T-^ f 



—A{v) , 

V 



r 



(50) 



U(l-a) 



V-41n(x) ,a = . 

We take the lower limit of integration in Xmin to be zero, since at low x (high v) 
the integrand is small. Then we have for R22 the following expression: 

R22 ~ 0A6GMsLc<7a^^ ■ (51) 

This conforms the value obtained in Ref. 7. Note that R22 does not depend on m 
and weakly depends on J, not depending on J at all in isothermal case. 
For Rii we have 



Ru 



GMsL,r+T-^i^-^ j J dxx^'^"i^{v) ■ (52) 



3V27r 



The integral diverges at lower limit, so we take the value (|41|l for Xmin- For the 
isothermal case (a = 0) it is possible to obtain an asymptotic approximation of 
this integral at low Xmin- Since r2o('S) ~ "7!^' integral approximately equals 
i—Xmin^'nXmin)~^ ~ 1//^- A similar consideration for a > leads to an estimate 
of the integral as The expression in brackets before the integral weakly 

depends on a and equals for a ~ 0. Then 



Rii « R22 ■ 0.1 i^-j . (53) 

The last coefficient R12 can be represented as R22 ■ Ka{fi), where Ka{fJ.) — > at 
fi ^ 0. So one can neglect the term with R12 in H25|l at low fi. 
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3.2. The diffusion coefficients for central region 

We restrict the calculation to the case when the density of stars n{r) cx r^-^/^. Then 




The integral weakly depends on Xmin at low Xmin and equals ~ 1.14. Finally, we 
have 



Comparing (|54|) with expression (|51|l for bulge, we conclude that they coincide at 
J = J(i?o) (i-e. at the boundary of central region of the black hole influence). This 
means that the expressions for two limiting cases are consistent with each other. 



Here the integral approximately equals 0.8/ ^Xmin', so with the aid of 147|) we obtain 

Ru - R22 ■ 0.25 - (56) 
m 

Similarly to the previous case we get that R12 — const R22 y , so it can be neglected 
at low m. 

3.3. Diffusion timescale estimate 

One can easily see that for an equation ^ — ^(^f^) the characteristic time of 

.2 

"equillibrium establishment" at spatial scale I equals t Let us make a few 

estimates for kinetic equation H17|) with coefficients (|51II53|) calculated for the bulge. 

The diffusion timescale for momenta T2 — ^g'^°^ ^ 10^ (l~pc) ^'^ ^^sX the 
galactic age corresponds to spatial area of r+ < 100 pc. The coefficient Ru increases 
with decreasing m. The diffustion timescale along J axis at minimal momentum 
m = nig equals ti ^ 3 • 10^ i^ic ^'^^ ^'^'^ ~ P*" ^1 much less that the 
galactic age. But for momenta at least an order of magnitude greater than nig the 
diffusion along J axis does not disturb much the one-dimensional diffusion along m. 
We will show in the next section that if we account for a cutoff of star distribution 
function at r ^ 0, then we obtain a finite limit for Ru with decreasing m. Thus 
Ti becomes comparable with T2 or greater for all r+ < 100 pc. Hence in the first 
approximation we can neglect the diffusion along J in comparison with the diffusion 
along m, at least for not very small values of m. 

In addition we determine the characteristic timescale for diffusion in the central 
region. Since R22 ^ = lo^ni^ the timescale does not depent on initial value of 



R22 



2AG^/^ Ms Mbh^^/^ Ro^^'^ ic ■ 



(54) 




(55) 



momentum and equals T2 = 10' 



yr X 
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3.4. Diffusion coefficients for an improved star distribution 

Indeed, in deriving relations (|33II38|I we have neglected the fact that number density 
of stars in the vicinity of black hole may be significantly less than that follows from 
general power-law profile due to stars capture or tidal disruption by black hole. From 
(|52|l we see that it is the stars in the pericenter of the orbit that affect the coefficient 
Rii most of all. Hence for particles with low momenta these are the innermost stars, 
which number we may have overestimated. Now let us try to evaluate corrections 
linked to this fact. 

Firstly, consider the Coulomb potential region. Observations show 1201 that the 
nearest star's orbital axis approximately equals 1000 AU = 3 • 10~^ pc. But from 
^ we see that the number of stars inside a sphere of radius 1000 AU should be 
about 15. 

For correction we adopt that stars distribution function vanishes for energies 
E < Ecr — — ^^7^ • Thus the density profile ((HJ transforms to 

3/2' 



n(r) =no [ — 



3/2 



1-(1-^ 

Tcr 



(57) 



Comparison with observations gives the value rcr ~ 5 • 10 pc ~ lO'^r, 



a- 



Now let us calculate modified coefficient Rn. Following H55|) we obtain the ex- 
pression (with integration from to 1 in this case): 



16 



Rn^-^G^Ms^L,For+\E\^/^T{ _ (58) 



3V2 



2Xc 



> Xcr > 1 



Xcr 

For r_|_ < 1 pc, = -y < 0.1, we have Xcr > Xmin, hence in the whole central region 
of Coulomb potential the lower limit of integration Rn is given by Xcr- This means 
that i?ii is in fact is much less (approximately ^\/xcr/Xmin times) than calculated 

from (innj- 

Secondly, the same formalism in the region of self-consistent (particularly, 
isothermal) potential should include the statement that distribution function of 
stars vanishes for sufficiently low energies. For isothermal star distribution this cor- 
responds to well-known solution for isothermal sphere with core (i.e. no central 
cusp), and the core radius is of the same order that rcr- Unfortunately, we can 
make no estimation for core radius from observations, since their spatial resolu- 
tion is not enough. We may adopt the same value as calculated for our Galaxy 
(rcr ~ 10^ -f- 10"* Tg). One can show that this assumption changes the lower limit 
of integration in 1)52(1 to value Xmin ^ max{r-,rcr)/'r+- Thus Rn increases with 
decreasing momentum up to m = (10^ 10'^)mg and then reaches constant limit. 
At the same time R22 is not affected by this cut since the integrand in H50|l is small 
at small x- 
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4. Dark matter absorption and growth of the black hole 
4.1. One- dimensional diffusion approximation 

It was shown in previous sections that coefBcient R22 does not depend on m in the 
bulge outside the central region of BH influence; coefRcient Rn ^ R22 {-^Y 1 e < li 
and R12 ^ i?22y • Since at the initial moment m ~ IqJ ^ J and our scope of interest 
lies in the domain of low momenta, we can take m <^ J and disregard the term with 
Ri2- Furthermore, from the same arguments it follows that ^ ^ at the initial 
moment, and since Rn ^ R22 at m = ^O'^i then Rudf/dJ <C R22df/dm at not 
very low values of m. Finally, /(m, J, t) = at to = to^, hence ^ =0, and 

m—mg 

Rii§j is limited at to ^ 0. Also notice that in fact Rn itself is limited at to ^ 0, 
as noted previously. 

To summarize, in the first approximation we leave in (|25|l only the first term 
and rewrite the kinetic equation ifTTjl as follows: 

J- = -R(J)—(m-L) , 59 

ot TO dm \ dm J 

d 

with boundary conditions f\m.=mn = : m——f\m=oo = , (60) 

om 

initial conditions 112(1 and diffusion coefRcient 

1 , J > Jo{Mbh) 



R{J)^0A6GMsLcao^^ ■ { f ^ ^ (61) 

J < Jn . 



Here Jq = ^GMth Rq/2 cx Ahh ^ (62) 

is the boundary value separating the BH domain of influence and bulge itself (as 
follows from [131)1. 

The flux of dark matter through the surface m — mg is given by the expression 

S{t) = -{2TTf III dJdmdm, % ^ -{^irf //" dJdm2m—R^m^ = 
.//./ ot II TO dm dm 



(27r)3 dJ -2 [ mR{J) 



^2{2nf dJhJ^'^Sj{t), (63) 



9/ 
dm 

where Sj{t) — mgR-^ is the flux through m = mg for equation H59|l with initial 
condition 

f{m, < = 0) = S{m^ - TOo^) , mn = loJ (64) 

Now we are going to calculate the time-dependent flux Sj{t) of DM particles 
through the absorption boundary for each value of J separately, i.e. one-dimensional 
flux, and then integrate it over J assuming that the diffusion along J is small. 
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4.2. Flux in one- dimensional diffusion 

Consider an auxilliary task: equation (|59|l with boundary conditions l|6()|l and initial 
condition (|64|) and determine the flux Sj{t) = nig R§^- 
The solution of ea. (|59|l may be represented as 

/"OO 

/(m,t)=/ dm' G{m,m',t) f{Tn',0) , (65) 
Jo 

G~ dXm' exp{~XRt) Z\{m,nig) Z\{m' ,mg) is the Green function , 
Jo 



Zx{m,mg) = 



iJo\VXmg) + Yo^{VXmg)y/^ 



is the orthogonal system of fundamental functions of the boundary problem H()()|> . 
Jo, 1^0 ^re Bessel functions of first and second kind of 0th order. 

The initial conditions give /(m, t) = 2:^G{ni, mo, t). 

One can easily show that 

d 



Q^Z,{m,mg) 



Then the flux 



Sj{t) = rUg R 



2mo 



T^rug (Jo2(v^mg) + V(VAmg))i/2 ' 

dZx{m,mg) 



dXnio exp{—XRt) Zx{mQ,m,g) 



dm 



RH{x,y) 
■nRt 



drj exp(— 77) 



Jo\^lx) + Yo\^x) 



mo N 
CRt' 



To obtain exact 



Here we have changed the variables: t] — XRt, x — mg/y/Rt, y — m^jy/ Rt. 

One can show that H{x,y) ~ Z{x,y) exp{—(^y'^) at y > x + 4, Z{x,y) weakly 
depends on its arguments, C ~ 5. Hence Sj{t) oc j exp(- 
form of the dependence, we take the following consideration. 

The flux Sj{t, m) — mR ^ is a continuous function of to; in the region mg < 
TO < Too we may take it to be a constant independent of m (it is correct for values 
of t greater than certain when the width of peak of DF becomes comparable with 

dfim,t) 



Too, see figEJ. Let us denote K{t) = 



dm 



, then we have 



Sj(t) ^ mg RK(t) , f{m,t) 



Rm 



■J dm' = mg K{t) In 



m 
m„ 



(66) 



Thus we obtain that f{m,t) grows logarithmically with to. Its dependence on 



time can be found from the argument that n{t) w jcxp^— ^^j, ^ is a con- 
stant, (fig. OJ- K{t) reaches maximal value at the exponent argument equal to — 1; 
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its maximal value K„iax — ^n^^ • '^'^ determine the value of S we apply expres- 
sion H66|) for m = mo/2; having found analytically fmaxifno/'^) = 2mo^ ' obtain 



e In 2 

Finally, the flux is 



The correctness of the above consideration is proved by numerical investigation 
of the problem. We also should notice that, in fact, expr. (I67II slightly overestimates 
the flux at t > (up to a factor of 1.5, see fig. Ej. 



4.3. Black hole growth law 

Firstly, we neglect that in central region of Coulomb potential the coefficient R{J) 
differs from the expression and take i?( J) = i?e J*^, e = < 0.2. Substituting 
the obtained value Hti7|l of Sj{t) into the expression Ht)3|l for total flux, we rewrite 
it as follows: 

S{t) = /°°/oJ^/«-^iexp ■ (68) 



In iai t V 5^^^ 

2mg ^ 

The absorption boundary nig changes with time. However, the flux weakly depends 
on the value nig. We approximate the logarithm in the denominator to have a 
constant value ~ 10. 



R 8(2-.) / 5 \ 8(2-.) 1 ( 



9 



Sit) = 8.9/og r:: . , = ( 7T ) —fcL . (69) 

Supposing that black hole growth is governed only by the absorption of dark 
matter, we obtain the following growth law for BH mass (assuming the mass of the 
seed black hole to be small): 

Mi,h{t)^B{R,t)^ , B = 8.9/oF,fci (70) 

Notice that more precise expression for Sj{t) at large t reduces the value of B 
approximately 1.2 times. 

Thus the black hole growth is power-law with power index about 9/16. This is 
in good agreement with previous work^l^ though the power index is a bit lower. 



4.4. Influence of central region onto the black hole growth 

One could suppose that the diffusion goes slower in the central region since the 
diffusion coefficient is lower there, and it may affect total growth law. In fact, 
nevertheless, it is not true. As has been shown in the characteristic timescale 
for diffusion r cx Mbh^ and is about 10^ yr. at the moment, hence it was even lower 
previously. Additionally, from 1)67(1 one can see that maximum of the flux goes from 



16 



International Journal of Modern Physics A 



inside the region Jmax 1^ y 5Rt/lo^, which at the moment corresponds to spatial 
area r < 100 pc, and value of Jq separating the central and outer parts of the bulge 
corresponds to r ~ 1 pc. These two quantities depend on time in similar ways: 
J7nax ^ t^^^, Jo ^ Mbhit) ^ t^/^^, so the relation Jq ^ Jmax was also true in the 
past. 

One can show that if R varies in time, then we should take the following expres- 
sion for Sj{t): 

Sj{t) = — ; — — exp j-^ 

For each J we define tmax{J) to be the time of maximal flux, and to (J) to be the 
time to enter the BH domain of influence: Jo(to) = J- The diffusion coefficient R{t) 
starts to decrease after the time to. But one can easily see that to{J) — 10"*^^^ = 
W^tmax{J), from which we conclude that this correction affects only the far "tail" 
of Sj{t), when then most part of the flux is already absorbed, and is practically 
unimportant for growth law. 



4.5. Influence of two-dimensional diffusion onto the black hole 
growth 

Now let's try to estimate the effect of diffusion along J axis, i.e. the correctness 
of one-dimensional approximation. We restrict our consideration to the case a = 
when R22 = const. Having obtained approximate solution for one-dimensional 
diffusion equation, we substitute it into initial equation H17|l and calculate the first 
term 



a / a/\ a / ^ J a [ e ln2 in- 



= aj r"ajj = aj l^^^^-^m aJ 1 1^1^!?^^^^ 



2 t2 1 



lo'J 



5R22t 



2m g 

Of special interest is the time t ^ tj^ax corresponding to the maximum of the flux: 
^max — \ji • It appears that t] oc t^f^ax 

51ni^ lo^Pm \ tmax J 

Now we should compare rj with the second term in H17|) rjm = -^-^ (^'n^R22§^^ ■ 
However, one can easily see that substituting H66() into this term makes it zero, 
since this approximate solution does not satisfy initial equation in the whole region 
TO > TOg. To avoid this problem, we calculate the left-hand side of (|59f) rjm = 

51n2 R22\n^ ^ 



T] 2 lo^ J 



Comparing with (|71(l we obtain — 



25 m 
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For J(r_|_ — 100 pc), m = mg this ratio is about 30. However, it was shown at 
the end of ii3.4l that to estimate maximal value of R\\ one should take m ~ lO^mg 
instead of m = m^, which lowers this ratio to about unity. 

It should be emphasized that all these estimations are rather approximate, since 
for precise assertions one should know the exact solution of one-dimensional diffu- 
sion. The arguments of H4.2I give only the correct value for flux through m — nig , 
but not the exact solution for all to. 

In general case, the inapplicability of reduction to one-dimensional diffusion at 
TO — > nig does not change much the value of flux through to = rUg, since the diffusion 
along J axis leads to the "blur" of distribution function along this axis, while the 
coefficient R22 only weakly depends on J as follows from (|53|l . This is a benefit of 
selecti ng r adial action as the second variable instead of energy as it was in previous 
works However, an especial consideration is necessary whether this diffusion 
may lead to the drift of dark matter into the region of large J where changes of 
TO, J during one period exceed their values, and diffusion approximation becomes 
incorrect. 

5. Comparison with observations and conclusions 

In conclusion, we make theoretical estimations of black hole masses obtained from 
growth law H70(l and compare them with observational data for several galaxies. 

As it was noted previously, in distant galaxies it is difficult to measure precisely 
the dependence of velocity dispersion on radius. However, in galaxies M 31 and 
NGC 4258 it seems to be almost constant and equals approximately 200 km/s^SHl 
Taking the value of /o the same as in our Galaxy H12(l and the time of growth 
t = 3 • 10^^ s = 10^" yr. we obtain the value Mbh = 1.8 • lO^Mg. The observed black 
hole masses are (2.0-^-8.5) • lO^M© for M 31 and 3.8 • IO^Mq for NGC 4258^. The 
comparison shows that dark matter may comprise significant fraction of black hole 
masses in these galaxies. 

As for Milky Way, the rotation curve is not flat in the center of the bulge, and 

the velocity dispersion may be approximately represented as cr(r) — ctq I -^^'p^ 1 , 



(To = 60 km/s^^. This corresponds to the value a = 0.5. Using the assessment (|12|l 
for /o the expression H70|l gives the black hole mass Mmi ~ IO^Mq, which clearly 
overestimates the adopted value of Mbh ~ 2.9 • IO^Mq -W' (about a factor of three). 
Notice that these values are about twice as smaller as in Ref. Q because of more 
precise estimate of dark matter flow H69|l . If we take a = const = 100 km/s, the 
numbers change only a little. 

The disagreement with observations may be explained by rather rough estimates 
of quantities /q and Iq, which are determined by the whole dark matter halo and 
are different for different galaxies. If we take the value for DM fraction in bulge 
T - 0.24 ^ 0.31 11^, then the resuhing BH mass will be almost in agreement with 
observation. 
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In conclusion, we can say that the model discussed can give reasonable estimate 
for observed masses of giant black holes in galactic centers. It is likely that a large 
fraction of black hole mass may be comprised of dark matter. 

Further development of this problem will require, firstly, more precise calculation 
of diffusion coefficients based on detailed data for star distribution in the central 
parts of the bulge; secondly, taking into account the bulge evolution; and finally, 
an exact consideration of two-dimensional diffusion accounting for particles Fermi- 
heating and their drift out of bulge. 

The authors are grateful to A.V. Gurevich, K.P. Zybin, A.S. Ilyin and 
V.A. Sirota for numerous fruitful discussions. This work was supported by RFBR 
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1 2 3 4 5 

V 



Fig. 1. The graphs of functions A{v) (coefficient A, EHJ)i ^i'") (coefficient A - B, gHJ) for 
different potentials: 1) power-law ('I' oc r"), a = 1/2; 2) power-law, a = 1/4; 3) isothermal 
(o = 0, * oc Inr). 
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Fig. 2. Normalized values f{m,t)/\n{m/mg) for different m 1661 : 1) m = r?io/2, 2) m = 20mg, 
3) m = 1.2mg. One can see that the graphs are quite similar, the normalized maximal values are 
approximately equal, though they are reached at different t. 
The values are given for mo = lOOrrig . 
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Fig. 3. The flux function K(t): 1) theoretical I67i . 2) obtained by numerical integration. One can 
see that the theoretical approximation is good enough for maximum of the flux, but overestimates 
it at large t. 
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Fig. 4. The corrected flux function K(t) for large t: 1) theoretical function (fig. EJ, 2) obtained 
by numerical integration, 3) correcter theoretical function; k oc t~^^^. The last approximation is 
much better. 
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